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We discuss the relationships between large deviations in stochastic systems, and “effective inter¬ 
actions” that induce particular rare events. We focus on the nature of these effective interactions 
in physical systems with many interacting degrees of freedom, which we illustrate by reviewing sev¬ 
eral recent studies. We describe the connections between effective interactions, large deviations at 
“level 2.5”, and the theory of optimal control. Finally, we discuss possible physical applications of 
variational results associated with those theories. 


I. INTRODUCTION 

Rare events are important in many physical settings: classic examples include phase transformation, protein-folding, 
and chemical reactions M- In those cases, a system makes a transition between two distinct states, and a variety of 
analytical and computational tools are available 0-Q. Here, we focus on a different class of rare events, where systems 
behave in an unusual fashion over an extended period of time. Specifically, we consider the probability of trajectories 
in which time-averaged quantities remain far om their typical (equilibrium) values. If the system is ergodic, the 
probabilities of such events decay to zero as the l eng th of trajectory goes to infinity: the rate of this decay is described 
by the mathematical theory of large deviations [9j. Recent studies of these larg e deviations have provided insights 
into fluctuation theorems [jjj, [lh], glassy systems [12l - fl5j , protein-folding [16H18I , chaotic dynamical systems [19j, I20| 
and interacting particle models |2l| - [24| . 

It turns out that the rare trajectories of interest in these systems can be characterised as typical trajectories for a 
certain modified system [25l - [30j , which we refer to here as the “auxiliary model”. The auxiliary model inherits many 
of its important properties from the original system of interest: if the original model has the Markov property then 
so does the auxiliary model. In many cases, the auxiliary model inherits the symmetries of the original model, and 
other properties like kinetic constraints are also preserved f20l . 

The existence of this auxiliary model raises important questions for characterisation of rare events. In particular, 
it means that by adding a particular set of interactions to the original model, one may drive the system to realise 
these rare events. In fact, these interactions can be shown to be the “optimal” ones for realising the rare events of 
interest, in a certain precise sense fjlffull (see Section HV Al below). It is therefore of great interest to characterise 
these interactions. For example, in glassy systems, they can stabilise “amorphous solid states” 00 that are 
otherwise only metastable - the nature of the interactions required to achieve this is a long-standing question in the 
field. In protein-folding systems, effective interactions might stabilise the native state, or they might favour misfolded 
states [13, m '■ understanding how these states can be characterised (and suppressed) is of vital importance in that 
context. 

In this paper, we survey some key results that are related to the existence and nature of these auxiliary models, 
and the effective interactions that they encode. Our aim is to draw together ideas from several different contexts and 
to give a (non-rigorous) presentation that highlights the central outstanding questions, and possible routes to solving 
them. In Section [H] we describe the setting for our main results. Section Hill illustrates the kinds of phenomena that 
we are interested in, through a summary of some recent numerical results. Then, in Section IIV1 we describe some 
theoretical results, including the relationship to large deviations at “level-2.5” j35| and to optimal control theory 0- 
134} . These results have not yet been exploited very far in the physics context - we highlight possibilities for future 
progress along these directions. Section fVH gives a brief summary and outlook. 


II. BASIC THEORY 

In this section, we collect some key results related to large deviations in stochastic processes. Many of these 
results have been derived independently in different contexts and by different groups. Here we follow the presentation 
of [H, [H, iHI; further details and references can be found in those works. 
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A. Models and master equations 

We consider a Markov process in continuous time, on a (finite) discrete state space with configurations {C}. For 
example, one can consider a lattice of Ising spins, or a simple particle model such as the asymmetric exclusion process. 
In practical settings, one is often interested in the thermodynamic limit, where the size of the state space is taken 
to infinity, for example by considering spins on lattices of increasing size. Alternatively, one may consider diffusive 
processes described by Langevin equations (stochastic differential equations). These may typically be obtained from 
lattice models by a continuum limit: one defines a process on a discrete lattice and then takes the lattice spacing to 
zero, rescaling time in an appropriate way to ensure diffusive behaviour. Our restriction to finite state spaces means 
that the following analysis may not always be valid on taking thermodynamic or continuum limits - in typical cases 
we expect our results to remain valid in such limits, but this is not guaranteed. 

The transition rates between configurations of the system are W(C' <— C). Let P(C,t ) be the probability that the 
system is in configuration C at time t: this quantity evolves by a Master equation 

d t P(C , t) = —r(C)P(C, t) + J2 W{C <- C')P{C', t ) (1) 

c 


where r(C) = Ylc W(C' <— C ) is the “escape rate” from configuration C. We assume that the process is irreducible, 
which ensures ergodicity, since the state space is finite. It is also useful to identify the subclass of these models that 
obey “detailed balance”. For these models, there exists a “potential” Ec such that 

W(C «- C')e~ Ec ' =W(C'<- C)e~ Ec , (2) 

for all C and C . Models with this property have time-reversal symmetric (“equilibrium”) steady states, in which the 
probability distribution over configurations is p(C) oc e~ Ec . 


B. Large deviations and biased ensembles 

The rare events that we consider are defined by the choice of an observable, which may be one of two types. 
A trajectory of the system consists of the (ordered) set of states which the system visits, and the times at which 
transitions (jumps) between states take place. The first type of observable takes the general form 

A = a (C' Q (3) 

jumps C—>C' 


where the sum runs over all transitions within the trajectory and the a(C’ C ) are a given set of numbers. For 
example, if a = 1 for all pairs of configurations then A is the total number of configuration changes in the trajectory. 
The second type of observable is the time integral of a state-dependent quantity 

obs 

B= dtb(C(t)). (4) 

Jo 

where C(t) is the configuration of the system at time t. For large f obs , the probability distribution of B generically 
has a large deviation form: 


p(B) ~ exp[-t obs </>(B/t obs )] (5) 

where cf>(b ) is known as a rate function. A similar expression holds for the distribution of A. [The precise meaning 
of (O is that lim toba _ >00 t~£ a h\p(B = bt 0bs ) = cj)(b); the symbol is used in this sense throughout this article.] The 
main question of interest in the following is: what kinds of dynamical trajectory dominate the distribution p(B) when 
B is not equal to its typical (steady-state) value? 

To obtain information about these trajectories, it is convenient to write a biased probability distribution over the 
possible trajectories of the model: 

- S B[C(t)] 

Prob[C(f); s] = Prob[C(f); 0] • ^ (6) 

where Prob[C(f);0] is the unbiased (steady-state) probability distribution over trajectories C(t), the parameter s sets 
the strength of the bias, and Z = (e~ sB )o resembles a partition function. Hence, the average any observable O within 
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this generalised ensemble is 


<o} s 


(Oe ~ sB ) 0 

Z (s,fobs) 


(7) 


It may be shown [§73] that averages within this biased ensemble are the same as those in an ensemble in which the 
value of A (or B) is constrained to a particular value. [Note that this equivalence is assured only in systems with 
finite state spaces, in which case the free energy i/j(s) is analytic and convex. In systems with infinite state spaces, 
dynamical phase transitions [TE [§6l . [52 ] may mean that trajectories which are representative of some values of A (or 
B) cannot be obtained within biased ensembles of the form given in (ftJ]) .] 

To analyse these biased ensembles, one considers the probability that a system is in configuration C at time f, and 
that the observable B has a particular value associated with the trajectory up to time t [ill . (MU- (The analysis 
for observables of type A is similar.) If this probability is p(C,B,t) then we define p(C,s,t) = J dBp(C, B, t)e~ sB . 
This quantity evolves by an equation which is formed from © by replacing P(C, t) with P(C, s, t ), and adding a term 
— sb{C)p(C , s, t) to the right hand side. The resulting equation is linear in P so it is useful to write it formally as 


d t \P) =W(s)|P) (8) 

where W (s) is an operator (matrix) with diagonal elements — r(C) — sb(C ) and off-diagonal elements W(C •<— C ). 
In the case of ty pe - A observables, the parameters a(C' C) appear in the off-diagonal elements via multiplicative 
factors e -s “ l§i,|3g|. Note that ([5]) resembles a master equation, but it does not conserve probability (in the sense 
that p(C, s, t) n °l constant under the time evolution). 


C. Connection between type-A and type-B observables 

We note at this point that the operator W(s) fully specifies the probability distribution in (0, up to possible 
boundary terms that we will neglect in the following (see also Sec. Ill Dl below). This means that if two processes are 
associated with the same operator W(s), then they have the same behaviour. It follows that ensembles defined by 
type-A observables can be given alternative definitions in terms of type -B observables, but for a different underlying 
stochastic model. 

For example, suppose that a model has transition rates W(C ■*— C ) and is biased by an observable of type A. Then, 
the same operator W(s) can be obtained by considering a different model with transition rates W{C <r- C) = W(C' 4— 
C)e -S “( c *~ c \ biased by an observable B = s -1 / d t[r(Ct) — r(C t )\ where f(C) = Yhc W(C' <— C): see for example [l3j, 
Appendix B]. A similar transformation means that any U-biased process can always be re-written as an A-biased 
one. (This requires that r(C ) + sb(C) > 0 for all configurations, which in finite state spaces can always be achieved by 
including an appropriate constant shift in b(C).) Hence, in the following, we sometimes state results either for type-A 
or type-B observables, since the results for the other type can always be derived by an appropriate transformation. 


D. Auxiliary models 

Given a model [specified by rates W(C' <— C)\ and an observable [specified by the a(C' <— C ) or 6(C)], one may always 
define an auxiliary model whose steady state distribution of trajectories is close to [A precise characterisation 
of this “closeness” is given in (flTl) below]. For observables of type B , the transition rates of the auxiliary model 

are ® miM® 

W am (C' <- C) = u C 'W{C <- Cjugi (9) 

where the uc are obtained by solving an eigenvalue equation for the operator W(s). Specifically, (u| is the left 
eigenvector associated with the smallest eigenvalue of — W(s): 

(it|(—W(s)) = i/)(s)(u\. (10) 

Here il>{s) is a dynamical free energy, related to the dynamical partition function by Z(s,tobs) ~ e~ tohB ^ l ' s ' > . We note 
the connection of 0 to Doob’s 6-transform ®]j which is one of the earliest results connecting rare events to auxiliary 
models of this kind. Similar results also appear in other kinds of biased rare-event problems [H, HH,®]- 
Eq. 0 motivates us to define an “effective potential” 


AVc = —2 lnuc. 


( 11 ) 
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With this definition, VF aux (C' <— C) = W(C' C)e^ AVc_Ayc '^ 2 , which can be interpreted as a modification of the 
original transition rates according to the change of the effective potential in a transition. For type-A observables, the 
analogue of © is 

W aux {C' ^C) = u C 'W(C «- C)e- sa{ - c '^ c) uc l (12) 

To see the relation between W aux and W(s), we define an operator W aux whose off-diagonal elements are the 
W aux (C' <— C ) and whose diagonal elements are —r aux (C), with escape rates r aux (C) = fF aux (C' C). If we also 

define u to be a diagonal operator whose elements are the uc, it follows [29| that 

W aux = nW(s)w -1 + ip, (13) 


which holds for both type-A and type-!? observables. 

With these definitions, the trajectory measure for the steady state of the auxiliary model, Prob[C(t)]; aux], is related 
to the biased ensemble © as 


e [ A Wt obs )-Ay C (o )]/ 2 

Prob[C(f);s] = Prob[C(f); aux]-- 


Po(C(0)) 
p aux (C( 0)) 


(14) 


where Z aux is a normalisation constant, and the final factor on the rhs is the ratio of steady state probabilities of the 
initial configuration C(0), in the original process [po(C(0))] and the auxiliary process [p aux (C(0))]. Eq. (fl4l) is most 
easily derived via direct construction of the various Prob[C(f)]. For example, if the biasing obervable is of type A then 
we have 


Prob[C(i), s] 


" K 

./:=i 


e -(to b s-tK)r(C K ) po ( Co j 


1 

Z(s, fobs) 


(15) 


where the trajectory is composed of configurations Co,Ci,... ,Ck , with configuration changes at times t\, £ 2 , • • • tic, 
we define to = 0, and po(C) is the probability of finding configuration C in the steady state of the original (unbiased) 
model. A similar construction of the analogous probability density for the auxiliary process then yields (1141) . For a 
detailed analysis, see [30]. 

Since all differences between the probability distributions in (fl^l) come from the initial and final states, we expect 
that for long trajectories, the distributions Prob[C(f);s] and Prob[C(f);aux] will differ only through initial and final 
“transient” regimes, and that their behaviour will be the same in the intermediate-time regime for which f > 1 and 
fobs - f > 1- The transient regimes are discussed in more detail in [Tt| | and also in [30|, where it was shown how a set 
of time-dependent auxiliary rates can lead to exact correspondence between the auxiliary and biased processes. 

It is useful to note that (for type-!? observables) 


u(C) oc lim (e sS + t ° b s’A s )) Cj0 (16) 

^obs ^00 

where the average is taken with respect to the unbiased dynamics, for a system initialised in configuration C [2(1 [27l . 
Hil, [58j]. The term t Q bs'i^(s) in the exponent ensures that the average does not grow or decay exponentially in time, 
because from the definition of Z and its link to the dynamical free energy one has 

e -tobsb(s) ^ Z( S: t ohs ) = {e~ sB )o- (17) 


E. Biased ensembles with time-reversal symmetry 


In cases where the biased ensembles are symmetric under time-reversal, the eigenvalue problem m may be sim¬ 
plified: it reduces to finding the largest eigenvalue of a symmetric matrix. The most common situation in which this 
occurs is when the unbiased model obeys detailed balance, and the biasing observable is either of type-!?, or of type-A 
with a(C' <— C) = a(C <— C) for all C and C. In this case one has simply EMI] 


ib = min 

\x) 


(s|e-®/ 2 (—W(s))e -®/ 2 |x) 

iA x ) 


(18) 


where E is a diagonal operator whose elements are the energies Ec that appear in the detailed balance relation ©, 
and the maximisation is over vectors with elements Xc■ The maximum occurs when xc = uce~ Ec / 2 so this variational 
result allows direct estimation of the effective interactions. Generalisations of this result to cases without time-reversal 
symmetry will be discussed in Section [IVI below. 
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III. ILLUSTRATIVE RESULTS FROM MODEL SYSTEMS 

Having introduced the general features of biased ensembles of trajectories, we now return to our original focus on 
complex systems with many interacting degrees of freedom. In these cases, it is not usually possible to solve the 
eigenproblem (1101) in order to obtain the uc■ Further, even if this eigenvector could be obtained exactly, it typically 
has such a large dimensionality that it does not provide direct information about the physical nature of effective 
interactions in the system. To illustrate these physical ideas, we now recall some recent results on the physical 
features of effective interactions in biased ensembles, for different model systems. 


A. Glass-forming systems 

Kinetically constrained models consist of interacting spins (or particles) in which local rules mean that only a subset 
of spins are able to flip at any given time step 0 , 0 These models provide simple descriptions of glass-forming 
liquids 0 ]. The “mobile” subset of spins changes with time, and the system is ergodic on long time scales. The 
dynamical motion in these systems can be complex and co-operative, even if their static (thermodynamic) properties 
are very simple. 

In these systems, it is typically possible to construct configurations in which the subset of mobile spins remains 
finite in the limit of large system size. In this case, if one considers the large deviations of the total number of spin 
flips in a trajectory (type- A observable with all a = 1), it can be shown from dl8l) that (i) limjv-^oo t/j/N < 0 where 
N is the system size, (ii) this bound is saturated for all s > 0, and (iii) the effective interaction in this case drives the 
system into configurations with a finite number of mobile spins. It follows that these systems have dynamical phase 
transitions at s = 0 [HI, [H|] . The dominant feature of the effective interactions for s > 0 is a very strong suppression 
of mobile spins, although the detailed nature of the effective interactions that produce this suppression is not known. 

Similar phase transitions exist in fully-connected (“mean-field”) spin-glass models with large numbers of metastable 
states 01 ] , and there is also numerical evidence for them in atomistic models of glass-forming liquids 00 , but the 
nature of the effective interactions again remains unclear. 

A recent study of a particular kinetically constrained model (the East model 0]) highlights the complex effective 
interactions that can appear even in simple systems. On biasing this model to low activity, one observes the dynamical 
phase transition discussed above. However, if one biases instead to high activity, one observes a hierarchy of responses 
that mirror the “aging” behaviour of the same model [0. (Aging behaviour occurs when the system is initialised at 
high temperature followed by dynamical relaxation at low temperature.) The dominant features of these states are 
(i) effective interactions that are long-ranged even for weak biases s, and (ii) a hierarchy of length scales associated 
with different relaxation processes within the system. 


B. Exclusion processes 

There have been many studies of large deviations in exclusion processes, in which particles move on a lattice, with 
at most one particle per site. Effective interactions in biased ensembles have been considered in relatively few cases; 
two examples are the limits of maximal dynamical activity or maximal current, where the effective interactions can be 
found exactly 0. These interactions are dominated by a long-ranged repulsion between particles: the system can be 
mapped to a “one-component plasma” of positively-changed particles interacting by Coulomb-like forces. The result 
of these long-ranged forces is that the system becomes “hyperuniform” 0] - density fluctuations on large length 
scales are strongly suppressed 0 ] . Such correlations occur in a variety of non-equilbrium systems 040 , but they 
are forbidden in equilibrium systems with short-ranged forces. 

Biasing exclusion processes to small activity can also result in phase transitions into inhomogeneous states 0, 
HH, H2|) although the effective interactions associated with these states have not been investigated in detail. Similar 
behavior can occur in simple models of heat conduction [0] . 


C. Numerical results 

As well as these analytic results, there are several numerical methods that allow large deviations to be investigated. 
Briefly, transition path sampling [7] is a computational method for sampling trajectories of systems according to 
general path ensembles, including examples such as © 00 - The method is most easily implemented for processes 
obeying detailed balance, although generalisations are possible [55| . Alternatively the cloning method was developed 
specifically to study large deviations 00,113 and is not restricted to systems with time-reversal symmetry - it 
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involves many copies (“clones”) of the system evolving in parallel. Finally, a third method was proposed recently 
by Nemoto and Sasa [58(, which involves direct estimation of the auxiliary rates in (1121) , in a manner reminiscent of 
thermodynamic integration. 

These methods have provided a number of interesting insights, especially for models that are not tractable analyt¬ 
ically. Examples include model protein-folding systems na, where biased ensembles are dominated by “misfolded” 
states, reminiscent of the low-activity states discussed in Sec. IIII A1 Similar results can also be obtained in protein 
systems for which Markovian effective descriptions are available - if the resulting state space is sufficiently small then 
large deviations can be analysed by exact diagonalisation of the operator W(s) |16j]- One again finds that the effective 
interactions stabilise misfolded metastable states fl8| . 

Numerical methods have also been used to study the competition between chaotic and periodic behaviour in 
dynamical systems mm. In particular, even if a system’s steady state is chaotic, its large deviations may be 
characterised by periodic trajectories, which allow the system to avoid “equilibration” into an ergodic state. 

We emphasise that the path sampling and cloning methods do not provide direct information about effective 
interactions, and even the method of [58j typically requires an approximate parameterisation of these interactions to be 
chosen before starting the analysis. However, the methods do yield representative configurations of the biased system, 
which at least provide qualitative insights into the underlying interactions. We believe that further development of 
methods in this area is a useful area for further study. 


D. General principles 

We identify two general principles from the illustrative examples above. Firstly, biased ensembles of trajectories 
often contain correlations that are very unusual in equilibrium systems. The liyperunifornr states found in exclusion 
processes are stabilised by long-ranged effective interactions [45|, |46| - these might not have been anticipated given 
the simple local rules and the simple bias to high activity. Similarly, the long-ranged correlated states found in the 
East model biased to high activity do not at all resemble the equilibrium state of that system [dij j. and nor do the 
periodic (non-chaotic) trajectories found in some dynamical systems [lj| [2(j. We emphasise that biased states are 
optimised with respect to global observables {A or B) that depend on the whole system, integrated over a long period 
of time, so there is no general reason to expect effective interactions to be the short-ranged forces that are familiar 
from equilibrium settings. So one may expect to find new and unusual phenomena on investigating large deviations. 

Secondly, effective interactions are often linked with underlying metastable states in a system biasing to low activity 
often drives the s yste m into “glassy” metastable states, as found in kinetically-constrained models [12], 13), atomistic 
glass-formers [l^llll: and proteins [l6l-[l8[. Given the variational principle (fT8[) . this may not be suprising - the low- 
lying eigenvalues of the operator — W(0) are naturally linked with metastable states and phase transitions, so weak 
perturbations can be expected to lead to hybridisation of these states with the dominant eigenvector. However, the 
use of large deviation methods to further analyse dynamical metastability and glassy behaviour seems promising. For 
example, recent work on biased ensembles in quantum systems also highlights the importance of quiescent (inactive) 
states that couple weakly to their enviroment {59J. 

IV. EFFECTIVE INTERACTIONS WITHOUT TIME-REVERSAL SYMMETRY 

This section surveys some results, mostly from the mathematical physics literature, which provide variational 
methods for determining uc in systems without time-reversal symmetry, so that (1181) does not apply. For systems 
of practical interest, we are proposing that these results could be useful for (i) analytic bounds on dynamical free 
energies (for example, proving the existence of phase transitions in non-equilibrium systems, following the analysis 
of the time-reversible case [l2|, HH); (ii) variational analyses of effective interactions, as used in (44} ; (iii) improved 
numerical procedures, for example obtaining an approximation to the auxiliary dynamics in order to improve sampling 
within a computational scheme. Our purpose here is to highlight these opportunities so we mostly quote relevant 
results, referring to the literature for more detailed analysis and derivations. 


A. Optimal control theory 

We first state a general variational formula for the free energy which may be viewed as a generalisation of 

(1181) for systems lacking time-reversal symmetry. The variation is over sets of transition rates, which should be chosen 
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to reproduce the auxiliary rates © as closely as possible. For type-P observables, 


V’(s) 


lim 

^ obs >-00 


1 

mm - 

{W™} t ohs 



L(C' 4 - C) 
VF var (C' <T- C) 


var. 


(19) 


where the variational parameters are (non-negative) rates W var (C / 4— C), the average is over a dynamical evolution 
under those rates starting from some arbitrary initial state, and 


L(C 4— C') = W vai (C 4— C') 


W var (C 4- C') 
n W(C<-C') 



+ W(C 4-C'). 


( 20 ) 


We note here an equivalent way of writing the objective function in (U9|) above. By averaging over the number of 
jumps in any small time interval after time t , starting from the current configuration C(f), one finds 


ip(s) 


lim 

^obs 


^obs 


mm - 

{W«“} fobs 


d 11 s b(C(t)) 


£L(C'<-C(t))\ 

V / v 


( 21 ) 


where we have also written out B explicitly as a time integral. The minima in (1191211) are obtained when the rates 
W var are equal to the auxiliary rates defined by ©. A derivation of this result will be sketched in Sec. IIV Bl below. 
We first give a brief discussion of its interpretation and potential usefulness. 

The variational principle (flT)l) arises in “optimal control theory” [3ll - f34| : the idea is that W var is a “controlled 
dynamics” that should be optimised in order to realise the rare event of interest. The content of m is that 
the controlled process should minimise s(B ) var , while deforming the original rates as little as possible. [Note that 
L(C 4— C ) resembles a relative entropy between the sets of transition rates, with L = 0 if W var (C' C) = W(C’ 4— C). 
In fact the final term in CD is exactly the small-At limit of the relative entropy between the distributions of 
configurations reached from C in a small time interval At when the rates W and VF var respectively are in force. For 
the rates W, this distribution is T\t{C r ) = A tW(C' C) for C ^ C and Pa*(C) = 1 — A tr(C) otherwise; the relevant 
expressions for the rates W var are analogous.] Since the maximum in (fTTJl) is obtained when W var = W aux , we may 
restrict the maximisation to rates W' ra *{C 4— C) = W(C' 4— C)eI AV var ( c )” A1 V “(C )]/2 0 f same f orm as w aux . Then 
AW ar has the interpretation of an effective potential that pushes the system towards the rare event of interest. In 
this context, (fl9l) can be interpreted as an optimisation over the “controlling field” AF var . 

In the case of diffusive processes, (1191) has a particularly simple form: consider a model defined by a Langevin 
equation (or stochastic differential equation) 


x = K(x) + r] (22) 

where K = K(x) is a force and rj is a white noise. We then define a “controlled process” x = K — d x V rvar + rj where 
V val = V val (x) is the controlling potential. The idea is to discretize in time using a small time interval At. For 
x' ~ x, one has 

TA/Warf ^ ^ ( x ' ~ x){~d x V) + exp((x' - x)d x V) - 1 (23) 

VV var (a:' 4— x) 

Then averaging over x' with weight W var reduces this to At(d x V) 2 /2 + 0(At 2 ). Hence 

4>(s)= lim min -^—f d t /sb{x{t)) + ^-[d x V{x(t))} 2 \ (24) 

toba^-oo y™ t obs J Q \ 2 / var 

where one seeks to simultaneously minimise the average of sB and the magnitude of the controlling force d x V var . The 
relationships between optimal control and large deviations for diffusive systems have been discussed in the physics 
literature [3j, |34[ , but while the results (1191211) for Markov chains are known in the mathematical literature [60|, they 
have not, to our knowledge, been applied very far in physics. 

In terms of future applications, it is clear that (1191) gives bounds on i/j and allows variational estimates of W aux . In 
principle this enables variational analyses of large deviations in non-equilibrium settings, similar to those described for 
time-reversible systems in Section fill Al However, there is an additional difficulty associated with (fl9l) . which arises 
from the estimation of the average with respect to the variational (controlled) dynamics. In the absence of detailed 
balance, these averages will typically need to be obtained by direct numerical simulation, in which case convergence 
to the limit of large t Q bs may be non-trivial. 



















In the case of time-reversal symmetric ensembles, one can restrict to W val that obey detailed balance, and © 
reduces to ©• To see this, replace the expectation value in (12T1) by an average with respect to the steady state 
of the controlled dynamics p val (C) oc q~ e ( c )~ av (O. The key point is that the logarithmic term in L{C' •<— C) 
yields £c,C' W™(C '<- C) p var (C)[AV {C) - AV(C')]/2; using the detailed balance relation W var (C' <- C)p vai (C) = 
W val (C <— C')p val (C') and interchanging the summation variables shows that this term vanishes. Finally using 
r(C) = £ c , W(C' ■<— C), Eq. (12T1) reduces to 


-0 


min ) 

{W /var } Z_^ 


sb(C) + r(C)~ E W var (C' <— C) 

C' 


M W (C) 


(25) 


which is the same as ©• 

We highlight two other potential routes for application of (1191) . First, it can provide simple bounds on ip by 
appropriate simple choices of W val . For example if one biases by the total activity (number of spin flips), and the 
system has a configuration with sub-extensive escape rate [there exists a sequence of configurations Cn in systems of 
increasing size N such that r(C/v)/-ZV —>■ 0 as N —>• 00 ], then linijv_> 00 ip(s)/N < 0 and hence (given weak conditions 
on properties of the steady state) there must be a dynamical phase transition at s = 0. This is a non-equilibrium 
analogue of results proven for kinetically constrained models of the glass transition [I 2 ], [l3j . It is relevant for exclusion 
processes, where the same result may be derived either by exact solution (52| or within fluctuating hydrodynamics [24} . 
But the method based on © is both very simple and very general. Second, there should be possibilities of using 
© in numerical schemes, for example by generalising the method of Nernoto and Sasa [58j . This possibility remains 
to be explored. 

Finally one could also consider finite-fobs analogues of m- We define <p(C,t 0 b s ) as the minimum value of the 
objective function on the r.h.s. of (1211) when starting from a given configuration C. It is then not difficult to argue 
that cp(C, fobs) = ip{s) + AVc/(2f 0 bs) for large fobs, up to corrections that decay exponentially with t 0 \, s . If one allows 
the variational rates W var to depend on time then one can also obtain a closed form for the evolution equation of the 
</>(C,f ob s). Thus it may be possible to obtain the effective interactions from the finite-fobs behaviour of the optimal 
cost <p(C,t 0 bs) in the control theory approach. 


B. Large deviations at “level-2.5” 

To understand the origin of the variational result ©, it is useful to consider the large deviations of a very 
general set of observables [28i] . For a given trajectory C(f), we define the empirical current which is a set of numbers 
Q(C <— C), obtained by counting the jumps (transitions) that between each pair of configurations, and dividing by 
fobs- Similarly, the empirical measure is a set of numbers p(C) given by the fraction of time that the trajectory spent 
in each configuration. Note that when we consider the number of jumps into and out of each configuration C along 
any trajectory, these two numbers must be equal or, if C is the initial or final configuration, differ by at most unity. 
Due to the division by fobs in the definition of Q(C <— C ) this difference vanishes for large fobs, leading to the balance 
condition 

^c) = J 2 q(c^c'). ( 26 ) 

a c 


1. Main result 


The observables p and Q are very high-dimensional objects if the state space is large, but as long as there are a 
finite number of them they obey a large deviation principle whose explicit rate function is known, for both biased 
and unbiased ensembles of trajectories. We first state the result for sufficiently large fobs, one has p(p,Q) ~ 
g — tobaA/bQ) with 


i{n,Q) 


E \ g( c ' <- c ) 

C,C' 


log 


Q{C <- C) 
W{C' <- C)p(C) 


- 1 


W{C' <- C)p{C) 


(27) 


In an ensemble biased by a type -B observable according to ©, all trajectories with a given p and Q are reweighted 
by the same factor e~ sB = exp (—s £ c b(C)p(C)), so after including the normalization factor l/e -tob8 ^( s ) one has 
p(p, Q) ~ e- t ohs[i(u,Q,s)-if>(s)] w Th 


I{p, Q , s) = I(p, Q) + s E KC)/r(C) 


( 28 ) 
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The result (177D is known as a “level-2.5” large deviation principle (LDP) since it is intermediate between an LDP 
for the empirical measure (known as level 2) and a full LDP for trajectories (known as level 3). A review of large 
deviations at level-2.5 is given in [35i |. while rigorous analysis of the case with countably infinite state spaces is given 
in [i61, [62], including a proof of (1271) . 


2. Connection to the auxiliary process 

The typical empirical measure and current within the biased ensemble of © can be obtained by minimisation 
over /i and Q of the rate function /(/it, Q , s) in (1271) . The key point for our purposes is that this allows a variational 
determination of the auxiliary model of (j9|) . We first cast the minimisation over p,, Q as a minimisation over a 
“variational auxiliary model”: a model for which typical trajectories have current Q and measure p,. The transition 
rates of this model then have to be 


W var (C' <- C) = Q(C '«- C)/n(C). (29) 

As one would expect, the balance constraint (1751) on Q then ensures that the steady state of the process defined by 
the rates W var is fi. 

We rewrite m as 


I(^Q) = J2\ wvai (C «-C) 

c,c 


j W var (C'<-C) ' 
° S W(C'<-C) 


+ W(C , «-C)j/i(C). 


(30) 


Including the bias term as in (1751) gives the function I(/i, Q, s), which is minimised by the (/i, Q) that are most likely 
within the biased ensemble: we denote their values by ( n*,Q*)■ The variational rates IU var at this minimum of 
I(fjL,Q,s) define a model for which typical trajectories have (/i, Q) = (/ P,Q *) so they must be exactly the auxiliary 
rates IU aux associated with the biased ensemble. 

Moreover, the minimal value of /(/r, Q, s) itself is just the dynamical free energy: ip(s) = I(n*, Q*,s ). This follows 
from a contraction principle 0 because the observable B = fobs He &(C)/i(C) * s a simple function of the empirical 
measure: recall from mu that (e sB )o ~ e Hence, decomposing the average into contributions from all 

possible one has 


( e ~ sB ) 


max 

n,Q 


e -t 0 bei(n,Q) e -stohs He b ( c )n(C) 


maxe _ 

n,Q 


(31) 


Summarising the ingredients so far, we have ip(s) = min Mi g J(/x, Q, s) where the minimization can equivalently be 
done over the variational rates W var rather than p, and Q. The final step in the argument is to realize that in the large 
f 0 bs - limit, the average over C{t) in the optimal control formulation (1211) above becomes an average over the stationary 
measure /i(C), making the penalty term L equal to I(n,Q) as rewritten in (130[) . Thus the variational principle (flTH) 
follows from the general level-2.5 result (1771) . 


3. Derivation of lw 


A derivation of (1301) for the case s = 0 is given in [28j . Here we outline their argument. The empirical current 
and measure (/it, Q) are typical for the process W var , but they are not typical for the original process W. For a given 
trajectory, one may obtain an expression for the ratio P[C(t); 0}/P[C(t); var], usng representations analogous to (fl5l) . 
Then one writes the probability of observing an empirical current and measure in the unbiased process as 


0 tobsl (MjQ) 


E p t c W;°] 

C(t)\n,Q 


E 


m*);0] 

P[C{t)\ var] 


P[C(t ); var] 


P[C(t); 0] \ 
P[C{t); var] / 


(32) 


The summations in this equation should be interpreted as path integrals over all trajectories that are compatible with 
an empirical current and measure (/it, Q). The average on the r.h.s. is with respect to the W var process: the restriction 
to a given (/it, Q) can be omitted here since this average is already dominated by such trajectories, which are typical 
for that process. Using the explicit form of the ratio to be averaged then yields (1771) . The analysis of [28| considered 
only the case s = 0, but the general result of (1751) follows immediately as explained above, because the effect of the 
bias in (Gj can be re-written as a bias that depends only on the empirical measure. 
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C. Large deviations at level-2 


Finally, we recall a classical result of Donsker and Varadhan 
Without bias, these satisfy p(p) ~ e -toh S J(n) with 


for large deviations of the empirical measure. 


JO) =maxj -J2 p(C')W(C' j 


^2r(C)/j(C). 


(33) 


where the maximisation is over a set of variational parameters p{C) > 0. In the B-biased ensemble the relevant large 
deviation function just needs to add the effect of the bias as before, giving p(p) ~ e _tobB w ^j 1 


J(p,s) = J(p) + ^sb(C)p(C). (34) 

e 

The corresponding expression for ensembles biased by type -A observables is given in [TH, Appendix C]. Subsequent 
minimisation over p yields the dynamical free energy ip, and the p(C) at the minimum are the u(C) associated with 
the auxiliary dynamics of Eq. ©• 

The result (l33l) can be obtained by minimisation of E71) over Q , subject to the balance constraints E6l) . Calling 
the minimum value JO), we want to show that it can be obtained alternatively from the maximisation problem (1331) . 
This can be done using Lagrangian duality: the Lagrangian for the original minimisation is 

LO, Q, A) = I(p, Q)+J2 A ( C )[Q( C ' <- C ) ~ Q( C <" <?')]■ (35) 

c,c 


The dual Lagrangian is then defined as LO, A) = rninQ L(p, Q, A). Since for any Q satisfying ECU) one has L(p, Q, A) = 
I(p,Q), it follows that L(p,X) < J(p). Since the equality holds for the optimal Q, one has the dual representation 
JO) = max a L O, A). Now setting the derivative of L{p, Q, A) to zero to find L(p, A) gives 


Q{C <- C) 
W{C' <- C)p(C) 


A(C') - A(C) 


(36) 


Substituting back into LO, Q , A), the log term cancels with the Lagrange multiplier contribution and one is left with 


LO, A) {~W(C t- C)p{C)e KC) ~ KC,) + W{C C)/i(C)} 
c,c 


(37) 


Identifying p(C) = e~ AI - c ) and carrying out the sum over C in the second term then gives (l33l) as desired. 

To our knowledge, Eq. (15H) has had limited application for estimation of and the u(C). One obstacle is that 
this requires a maximisation over p, followed by a minimisation over p. For this reason, straightforward bounds on ip 
are not directly available, unlike the case of (E71) where one minimises over both p and Q. 


V. OUTLOOK 

We have summarised a range of analytical and numerical results related to the effective potentials encoded by 0111- 
Section EH reviews some previous results where these effective potentials have been estimated, mostly in time-reversal 
symmetric ensembles. Section m shows how the effective potentials can be interpreted in terms of the controlling 
forces that achieve rare events most efficiently, in the sense of the “objective function” L in (EU1) . We have discussed 
how the variational results described in Sections IIV Al and IIV Bl might be useful for generalising these kinds of method 
to systems without detailed balance, and for developing new numerical methods, possibly following Ref. [58j |. The 
application of these results to biased ensembles for open quantum systems [59j might also provide useful insights. 

Another general challenge coming from biased ensembles is the description of biased states that are inhomogeneous 
in space and time. The “addivity principle” leads to some exact results in homogeneous systems, but an accurate 
description of spatially inhomogeneous (phase-separated) states remains outstanding in some cases. Biased ensembles 
also support “travelling-wave” states which are inhomogeneous in both space and time [52, [(33j: it might be useful to 
investigate variational techniques based on ED in order to address these problems. 

From a fundamental point of view, the relation between effective interactions and the thermodynamic limit is 
also important. Biased ensembles in general will be characterised by some stationary measure p(C). Restricting for 
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convenience to systems with time-reversal symmetry one then expects that this has the form ^i(C ) ex no(C)e~ AVc , 
where /i 0 is the stationary distribution of the unbiased process, and A Vc an effective potential. However, in the 
thermoydynamic limit, a question arises as to whether the measure fi is “Gibbsian” I64l l65|: that is, whether AH can 
be written as a well-defined sum of interaction terms of increasing range. If such a description is not possible, even the 
definition of effective interactions becomes problematic in the thermodynamic limit. If one considers large deviations 
of the total energy (type-73) in the Ising model, there is evidence that the resulting effective interactions may not be 
Gibbsian [29i | . It is also not clear how whether the limits of large system size and large-t 0 bs should commute in such 
cases, and what consequences this might have. It would be interesting to analyse these questions further in future 
work. 
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